import os
import getopt
import sys
import string
from parse_fasta_file import get_tag_seq_from_fasta_w_substr
from parse_fasta_file import reverseCompliment
'''
This scripts converts *.junctions file to *junc filr for TopHat
*.junctions has format: chr start end gene support strand.
*.juncs has format: chr start end gene support strand. 
'''
def junctions2juncs( iFileN, oFileN ):
    iFile = open( iFileN, "r" )
    oFile = open( oFileN, "w" )
    EXP_INPUT_COL_N = 6
    # SELECTED_OUT_COL = [0, 1, 2, 5]
    while True:
        line = iFile.readline()
        if not line:
            break
        info = line.split( '\t' )
        if ( len( info ) == EXP_INPUT_COL_N ):
            if len( info[2] ) <= 0 :
                continue
            exonStart = int( info[1] ) - 1
            exonEnd = int( info[2] )
            oLine = info[0] + '\t' + \
                    str( exonStart ) + '\t' + \
                    str( exonEnd ) + '\t' + \
                    info[5]
            oFile.write( oLine )
    iFile.close()
    oFile.close()

def checkJuncs( juncFN ):
    juncF = open( juncFN, "r" )
    EXP_INPUT_COL_N = 4
    NO_CHECK_JUNCS = 100
    [seqN, seq] = ["", ""]
    for i in range( NO_CHECK_JUNCS ):
        line = juncF.readline()
        if not line:
            break
        info = line.split( '\t' )
        if len( info ) != EXP_INPUT_COL_N:
            print info
            continue
        chrN = info[0]
        start = int( info[1] )
        end = int( info[2] )
        strand = info[3].strip()
        PRINT_JUNC_RANGE = 5
        if i == 0:
            iFile = open( chrN + '.fa', "r" )
            [seqN, seq] = get_tag_seq_from_fasta_w_substr( iFile, '>' )
            iFile.close()
            seqL = len( seq )
            print seqN, seqL, '\n'
        if ( len( info ) != EXP_INPUT_COL_N ):
            print "unexpected line\n", info
            continue
        if ( seqL <= end ):
            print "seq is too short\n", seqL, end
            continue

        JuncEISeq = seq[start - PRINT_JUNC_RANGE + 1: start + 1].upper()\
                + seq[start + 1: start + PRINT_JUNC_RANGE + 1].lower()
        JuncIESeq = seq[end - PRINT_JUNC_RANGE: end].lower()\
                + seq[end: end + PRINT_JUNC_RANGE].upper()
        if ( strand == '-' ):
            ( JuncEISeq, JuncIESeq ) = ( JuncIESeq, JuncEISeq ) # swap
            JuncIESeq = string.join( reverseCompliment( JuncIESeq ), "" )
            JuncEISeq = string.join( reverseCompliment( JuncEISeq ), "" )
        print JuncEISeq, JuncIESeq, i, strand

    iFile.close()
    juncF.close()

def usage():
    sys.stderr.write( "Working dir is " + os.getcwd() + '\n' )
    print __doc__
    print "usage: ", sys.argv[0]
    print "The 1st ptrs is the *.junctions file name."
    print "The 2nd ptrs is the *.juncs file name."
    sys.exit( 2 )

def main():
    try:
        opts, args = getopt.getopt( sys.argv[1:], 'h', ["help"] )
    except getopt.error:
        usage()
        sys.exit( 2 );
    # process options and argument 
    for o in opts:
        if o in ( "-h", "--help" ):
            usage()
            sys.exit();
    argc = len( args )
    if( argc == 2 ):
        junctions2juncs( args[0], args[1] )
        checkJuncs( args[1] )
    else:
        usage()
        print argc, "args are given"

if __name__ == '__main__':
    main()
